* ADL DGP 

cscript
* Set working directory
cd "/Users/alikagalwala/Library/CloudStorage/Dropbox/Whitten_Kagalwala/TS/JOP/Acceptance/KW_Replication/ADL_DGP"

* Run the file with the simulation program
do "ADL-DGP.do"

// Monte Carlo Analyses
* Create a matrix to store results
matrix mc_summary_stats = J(175, 71, .)

matrix colnames mc_summary_stats = t phi_y phi_x ///
phi_adl phi_adl_t1 phi_adl_power b1_adl b1_adl_t1 b1_adl_power b2_adl b2_adl_t1 b2_adl_power lrm_adl lrm_adl_t1 lrm_adl_power ///
phi_gecm phi_gecm_t1 phi_gecm_power b1_gecm b1_gecm_t1 b1_gecm_power b2_gecm b2_gecm_t1 b2_gecm_power lrm_gecm lrm_gecm_t1 lrm_gecm_power ///
phi1_adl22 phi1_adl22_t1 phi1_adl22_power phi2_adl22 phi2_adl22_t1 b1_adl22 b1_adl22_t1 b1_adl22_power b2_adl22 b2_adl22_t1 b2_adl22_power ///
b3_adl22 b3_adl22_t1 lrm_adl22 lrm_adl22_t1 lrm_adl22_power ///
phi1_adl33 phi1_adl33_t1 phi1_adl33_power phi2_adl33 phi2_adl33_t1 phi3_adl33 phi3_adl33_t1 b1_adl33 b1_adl33_t1 b1_adl33_power b2_adl33 b2_adl33_t1 b2_adl33_power b3_adl33 b3_adl33_t1 b4_adl33 b4_adl33_t1 lrm_adl33 lrm_adl33_t1 lrm_adl33_power ///
ftest_y_adl22_t1 ftest_y_adl22_power ftest_x_adl22_t1 ftest_x_adl22_power ///
ftest_y_adl33_t1 ftest_y_adl33_power ftest_x_adl33_t1 ftest_x_adl33_power

* Create a counter to store results in the matrix
local j = 1

foreach t in 25 50 75 100 150 500 1000 {
	foreach rho_e in 0 0.2 0.4 0.6 0.8 {
		foreach rho_x in 0 0.2 0.4 0.6 0.8 {
			parallel sim, expr(phi_adl = r(phi_adl)  phi_adl_t1 = r(phi_adl_t1) phi_adl_power = r(phi_adl_power) ///
			b1_adl = r(b1_adl) b1_adl_t1 = r(b1_adl_t1) b1_adl_power = r(b1_adl_power) b2_adl = r(b2_adl) b2_adl_t1 = r(b2_adl_t1) ///
			b2_adl_power = r(b2_adl_power) lrm_adl = r(lrm_adl) lrm_adl_t1 = r(lrm_adl_t1) lrm_adl_power = r(lrm_adl_power) ///
			phi_gecm = r(phi_gecm) phi_gecm_t1 = r(phi_gecm_t1) phi_gecm_power= r(phi_gecm_power) b1_gecm = r(b1_gecm) ///
			b1_gecm_t1 = r(b1_gecm_t1) b1_gecm_power= r(b1_gecm_power) b2_gecm = r(b2_gecm) b2_gecm_t1 = r(b2_gecm_t1) ///
			b2_gecm_power = r(b2_gecm_power) lrm_gecm = r(lrm_gecm) lrm_gecm_t1 = r(lrm_gecm_t1) lrm_gecm_power = r(lrm_gecm_power) ///
			phi1_adl22 = r(phi1_adl22) phi1_adl22_t1 = r(phi1_adl22_t1) phi1_adl22_power = r(phi1_adl22_power) phi2_adl22 = r(phi2_adl22) ///
			phi2_adl22_t1 = r(phi2_adl22_t1) b1_adl22 = r(b1_adl22) b1_adl22_t1 = r(b1_adl22_t1) b1_adl22_power = r(b1_adl22_power) ///
			b2_adl22 = r(b2_adl22) b2_adl22_t1 = r(b2_adl22_t1) b2_adl22_power = r(b2_adl22_power) b3_adl22 = r(b3_adl22) ///
			b3_adl22_t1 = r(b3_adl22_t1) lrm_adl22 = r(lrm_adl22) lrm_adl22_t1 = r(lrm_adl22_t1) lrm_adl22_power = r(lrm_adl22_power) ///
			phi1_adl33 = r(phi1_adl33) phi1_adl33_t1 = r(phi1_adl33_t1) phi1_adl33_power = r(phi1_adl33_power) phi2_adl33 = r(phi2_adl33) ///
			phi2_adl33_t1 = r(phi2_adl33_t1) phi3_adl33 = r(phi3_adl33) phi3_adl33_t1 = r(phi3_adl33_t1) b1_adl33 = r(b1_adl33) ///
			b1_adl33_t1 = r(b1_adl33_t1) b1_adl33_power = r(b1_adl33_power) b2_adl33 = r(b2_adl33) b2_adl33_t1 = r(b2_adl33_t1) ///
			b2_adl33_power = r(b2_adl33_power) b3_adl33 = r(b3_adl33) b3_adl33_t1 = r(b3_adl33_t1) b4_adl33 = r(b4_adl33) ///
			b4_adl33_t1 = r(b4_adl33_t1) lrm_adl33 = r(lrm_adl33) lrm_adl33_t1 = r(lrm_adl33_t1) lrm_adl33_power = r(lrm_adl33_power) ftest_y_adl22_t1 = r(ftest_y_adl22_t1) ftest_y_adl22_power = r(ftest_y_adl22_power) ///
					ftest_x_adl22_t1 = r(ftest_x_adl22_t1) ftest_x_adl22_power = r(ftest_x_adl22_power) ftest_y_adl33_t1 = r(ftest_y_adl33_t1) ftest_y_adl33_power = r(ftest_y_adl33_power) ftest_x_adl33_t1 = r(ftest_x_adl33_t1)  ///
					ftest_x_adl33_power = r(ftest_x_adl33_power)) ///
			reps(1000) noisily seeds(12344 12344 12344 12344 12344 12344 12344 12344): adl, obs(`t') rho_e(`rho_e') rho_x(`rho_x')
			  
			matrix mc_summary_stats[`j', 1] = `t'
			matrix mc_summary_stats[`j', 2] = `rho_e'
			matrix mc_summary_stats[`j', 3] = `rho_x'
			
			* ADL (1,1)
			qui su phi_adl
			matrix mc_summary_stats[`j', 4] = r(mean)
			qui su phi_adl_t1
			matrix mc_summary_stats[`j', 5] = r(mean)
			qui su phi_adl_power
			matrix mc_summary_stats[`j', 6] = r(mean)
			qui su b1_adl
			matrix mc_summary_stats[`j', 7] = r(mean)
			qui su b1_adl_t1
			matrix mc_summary_stats[`j', 8] = r(mean)
			qui su b1_adl_power
			matrix mc_summary_stats[`j', 9] = r(mean)
			qui su b2_adl
			matrix mc_summary_stats[`j', 10] = r(mean)
			qui su b2_adl_t1
			matrix mc_summary_stats[`j', 11] = r(mean)
			qui su b2_adl_power
			matrix mc_summary_stats[`j', 12] = r(mean)
			qui su lrm_adl
			matrix mc_summary_stats[`j', 13] = r(mean)
			qui su lrm_adl_t1
			matrix mc_summary_stats[`j', 14] = r(mean)
			qui su lrm_adl_power
			matrix mc_summary_stats[`j', 15] = r(mean)
			
			* GECM
			qui su phi_gecm
			matrix mc_summary_stats[`j', 16] = r(mean)
			qui su phi_gecm_t1
			matrix mc_summary_stats[`j', 17] = r(mean)
			qui su phi_gecm_power
			matrix mc_summary_stats[`j', 18] = r(mean)
			qui su b1_gecm
			matrix mc_summary_stats[`j', 19] = r(mean)
			qui su b1_gecm_t1
			matrix mc_summary_stats[`j', 20] = r(mean)
			qui su b1_gecm_power
			matrix mc_summary_stats[`j', 21] = r(mean)
			qui su b2_gecm
			matrix mc_summary_stats[`j', 22] = r(mean)
			qui su b2_gecm_t1
			matrix mc_summary_stats[`j', 23] = r(mean)
			qui su b2_gecm_power
			matrix mc_summary_stats[`j', 24] = r(mean)
			qui su lrm_gecm
			matrix mc_summary_stats[`j', 25] = r(mean)
			qui su lrm_gecm_t1
			matrix mc_summary_stats[`j', 26] = r(mean)
			qui su lrm_gecm_power
			matrix mc_summary_stats[`j', 27] = r(mean)
		
			* ADL (2,2)
			qui su phi1_adl22
			matrix mc_summary_stats[`j', 28] = r(mean)
			qui su phi1_adl22_t1
			matrix mc_summary_stats[`j', 29] = r(mean)
			qui su phi1_adl22_power
			matrix mc_summary_stats[`j', 30] = r(mean)
			qui su phi2_adl22
			matrix mc_summary_stats[`j', 31] = r(mean)
			qui su phi2_adl22_t1
			matrix mc_summary_stats[`j', 32] = r(mean)
			qui su b1_adl22
			matrix mc_summary_stats[`j', 33] = r(mean)
			qui su b1_adl22_t1
			matrix mc_summary_stats[`j', 34] = r(mean)
			qui su b1_adl22_power
			matrix mc_summary_stats[`j', 35] = r(mean)
			qui su b2_adl22
			matrix mc_summary_stats[`j', 36] = r(mean)
			qui su b2_adl22_t1
			matrix mc_summary_stats[`j', 37] = r(mean)
			qui su b2_adl22_power
			matrix mc_summary_stats[`j', 38] = r(mean)
			qui su b3_adl22
			matrix mc_summary_stats[`j', 39] = r(mean)
			qui su b3_adl22_t1
			matrix mc_summary_stats[`j', 40] = r(mean)
			qui su lrm_adl22
			matrix mc_summary_stats[`j', 41] = r(mean)
			qui su lrm_adl22_t1
			matrix mc_summary_stats[`j', 42] = r(mean)
			qui su lrm_adl22_power
			matrix mc_summary_stats[`j', 43] = r(mean)
	
			* ADL (3,3)
			qui su phi1_adl33
			matrix mc_summary_stats[`j', 44] = r(mean)
			qui su phi1_adl33_t1
			matrix mc_summary_stats[`j', 45] = r(mean)
			qui su phi1_adl33_power
			matrix mc_summary_stats[`j', 46] = r(mean)
			qui su phi2_adl33
			matrix mc_summary_stats[`j', 47] = r(mean)
			qui su phi2_adl33_t1
			matrix mc_summary_stats[`j', 48] = r(mean)
			qui su phi3_adl33
			matrix mc_summary_stats[`j', 49] = r(mean)
			qui su phi3_adl33_t1
			matrix mc_summary_stats[`j', 50] = r(mean)
			qui su b1_adl33
			matrix mc_summary_stats[`j', 51] = r(mean)
			qui su b1_adl33_t1
			matrix mc_summary_stats[`j', 52] = r(mean)
			qui su b1_adl33_power
			matrix mc_summary_stats[`j', 53] = r(mean)
			qui su b2_adl33
			matrix mc_summary_stats[`j', 54] = r(mean)
			qui su b2_adl33_t1
			matrix mc_summary_stats[`j', 55] = r(mean)
			qui su b2_adl33_power
			matrix mc_summary_stats[`j', 56] = r(mean)
			qui su b3_adl33
			matrix mc_summary_stats[`j', 57] = r(mean)
			qui su b3_adl33_t1
			matrix mc_summary_stats[`j', 58] = r(mean)
			qui su b4_adl33
			matrix mc_summary_stats[`j', 59] = r(mean)
			qui su b4_adl33_t1
			matrix mc_summary_stats[`j', 60] = r(mean)
			qui su lrm_adl33
			matrix mc_summary_stats[`j', 61] = r(mean)
			qui su lrm_adl33_t1
			matrix mc_summary_stats[`j', 62] = r(mean)
			qui su lrm_adl33_power
			matrix mc_summary_stats[`j', 63] = r(mean)
			
			* F-tests
			qui su ftest_y_adl22_t1
			matrix mc_summary_stats[`j', 64] = r(mean)
			qui su ftest_y_adl22_power
			matrix mc_summary_stats[`j', 65] = r(mean)
			qui su ftest_x_adl22_t1
			matrix mc_summary_stats[`j', 66] = r(mean)
			qui su ftest_x_adl22_power
			matrix mc_summary_stats[`j', 67] = r(mean)
			qui su ftest_y_adl33_t1
			matrix mc_summary_stats[`j', 68] = r(mean)
			qui su ftest_y_adl33_power
			matrix mc_summary_stats[`j', 69] = r(mean)
			qui su ftest_x_adl33_t1
			matrix mc_summary_stats[`j', 70] = r(mean)
			qui su ftest_x_adl33_power
			matrix mc_summary_stats[`j', 71] = r(mean)
			
			* Iterate the counter
			local j = `j' + 1
		}
	}
}

clear
set obs 175
svmat double mc_summary_stats, names(col)
save adl.dta, replace

* Bias calculations
gen bias_adl_phi = phi_adl - phi_y
gen bias_adl_b1 = b1_adl - 2
gen bias_adl_b2 = b2_adl - 1.5
gen bias_adl_lrm = lrm_adl - ((2 + 1.5)/(1-phi_y))

gen bias_gecm_phi = phi_gecm - (phi_y - 1)
gen bias_gecm_b1 = b1_gecm - 2
gen bias_gecm_b2 = b2_gecm - (2 + 1.5)
gen bias_gecm_lrm = lrm_gecm - ((2 + 1.5)/(1 - phi_y))

gen bias_adl22_phi = phi1_adl22 - phi_y
gen bias_adl22_b1 = b1_adl22 - 2
gen bias_adl22_b2 = b2_adl22 - 1.5
gen bias_adl22_lrm = lrm_adl22 - ((2 + 1.5)/(1-phi_y))
gen bias_adl22_phi2 = phi2_adl22 - 0
gen bias_adl22_b3 = b3_adl22 - 0

gen bias_adl33_phi = phi1_adl33 - phi_y
gen bias_adl33_b1 = b1_adl33 - 2
gen bias_adl33_b2 = b2_adl33 - 1.5
gen bias_adl33_lrm = lrm_adl33 - ((2 + 1.5)/(1-phi_y))
gen bias_adl33_phi2 = phi2_adl33 - 0
gen bias_adl33_phi3 = phi3_adl33 - 0
gen bias_adl33_b3 = b3_adl33 - 0
gen bias_adl33_b4 = b4_adl33 - 0

save adl.dta, replace
